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Abstract 

We investigate a single particle on a 3-dimensional, cubic lattice with a random on-site potential (3D Anderson model). We 
concretely address the question whether or not the dynamics of the particle is in full accord with the diffusion equation. Our 
approach is based on the time-convolutionless (TCL) projection operator technique and allows for a detailed investigation 
of this question at high temperatures. It turns out that diffusive dynamics is to be expected for a rather short range of 
wavelengths, even if the amount of disorder is tuned to maximize this range. Our results are partially counterchecked by the 
numerical solution of the full time-dependent Schrodinger equation. 
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Since it had been suggested by P. W. Anderson, the 
Anderson model served as a paradigm for transport in 
disordered systems [1,2,3,4,5,6]. In its probably sim- 
plest form the Hamiltonian may be written as 



H = y e T a\ a r + > ^ a\ a r 



(1) 



where a r , at are the usual annihilation, respectively 
creation operators; r labels the sites of a d- dimensional 
lattice; and NN indicates a sum over nearest neighbors. 
The e r are independent random numbers, e.g., Gaus- 
sian distributed numbers with mean (e r ) = and vari- 
ance (e r e r '} = 5 T T i a 2 . Thus, the first sum in Eq. (1) 
describes a random on-site potential and hence disor- 
der. 

The phenomenon of localization, including localiza- 
tion lengths, has intensively been studied in this sys- 
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tem [1,2,3,4]. For the lower dimensional cases d — 1 
and d — 2 (in the thermodynamic limit, i.e., with re- 
spect to the infinite length scale) an insulator results 
for arbitrary (non-zero) values of a, see, e.g., Ref. [4] . Of 
particular interest is the 3-dimensional case. Here, at 
zero temperature, increasing disorder induces a metal- 
to- insulator transition at the infinite length scale [3,4]. 
However, with respect to finite length scales the follow- 
ing transport types are generally expected: i.) ballistic 
on a scale below some, say, mean free path; ii.) possi- 
bly diffusive on a scale above this mean free path but 
below the localization length; iii.) localized (isolating) 
on a scale above the localization length. In the above 
transition decreasing disorder is viewed to reduce the 
size of the third regime, until it vanishes. 

Here, other than most of the pertinent literature, 
we do not focus on the mere existence of a finite lo- 
calization length. We rather concentrate on the size 
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of the intermediate regime and the dynamics within. 
We demonstrate that it is indeed diffusive (rather than 
subdiffusive, superdiffusive, or anything else). In prin- 
ciple, for long localization lengths (or no localization) 
this regime could be very large. But the results pre- 
sented in the paper at hand indicate that it is not, at 
least not in the limit of high temperatures. Investiga- 
tions in that direction (but not for d — 3) are also per- 
formed in Refs. [6,7]. 

Our approach is based on the time-convolutionless 
(TCL) projection operator technique [8] which has al- 
ready been applied to the transport properties of simi- 
lar models without disorder [9,10,11,12]. In its simplest 
form (which we apply here) this method is restricted 
to the infinite temperature limit. This implies that en- 
ergy dependences are not resolved, i.e., our results are 
to be interpreted as results on an overall behavior of 
all energy regimes. 
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Fig. 1. A 3-dimcnsional lattice which consists of N layers 
with n X n sites each. Only next-neighbor hoppings are 
taken into account. Intra-layer hoppings are specified by 
a constant a = 1 (white arrows), inter-layer hoppings by 
another constant fi (black arrows). 

As shown in Fig. 1, we consider a 3-dimensional lat- 
tice consisting of N layers with n x n sites each. The 
Hamiltonian of our model is almost identical to Eq. (1) 
with one single exception: All hopping terms that cor- 
respond to hoppings between layers (black arrows in 
Fig. 1) are multiplied by some constant (5. This is ba- 
sically done due to technical reasons, see below. How- 
ever, for ft = 1 the Hamiltonian reduces to the stan- 
dard Anderson Hamiltonian (1). 

We now establish a "coarse-grained" description in 
terms of subunits: At first we take all those terms of 
the Hamiltonian which only contain the sites of the 
/ith layer in order to form the local Hamiltonian of 
the subunit p. Thereafter all those terms which contain 
the sites of adjacent layers \x and fj, + 1 are taken in 
order to form the interaction (5 ty between neighboring 



subunits fi and //+ 1. Then the total Hamiltonian may 
be also written as H = Ho + f3 V, 

N-l N-l 

Ho=J2 h »> ^ =$>#•■ (2) 

where we employ periodic boundary conditions, e.g., 
we identify p = N with fi — 0. The above introduction 
of the additional parameter f3 thus allows for the inde- 
pendent adjustment of the "interaction strength" . We 
are going to work in the interaction picture. The hence 
indispensable eigenbasis of Ho may be found from the 
diagonalization of disconnected layers. 

By n^, we denote the particle number operator of 
the /xth subunit, i.e., the sum of a].a r over all r of 
the /uth layer. Since [J^ n M , H] = 0, the one-particle 
subspace may be analyzed separately, which will be 
done throughout this work. 

The total state of the system is naturally repre- 
sented by a time-dependent density matrix p(t). Con- 
sequently, the quantity P M (i) = Tr{ p(t) n M } is the 
probability for locating the particle somewhere within 
the /ith subunit. The consideration of these "coarse- 
grained" probabilities corresponds to the investigation 
of transport along the direction which is perpendicular 
to the layers. Instead of simply characterizing whether 
or not there is transport at all, we analyze the full dy- 
namics of the P M . 

Those dynamics may be called diffusive, if the P M 
fulfill a discrete diffusion equation 



Pu = K(P„+i + P l ,-i-2P fi ) 



(3) 



with some diffusion constant k. A decoupled form of 
this equation is routinely derived by a transformation 
onto, e.g., cosine-shaped Fourier modes, that is, Eq. (3) 
yields 

N-l 

F q = -2 (1 - cos q) k F q , F q = C q cos(g fi) P M 

fi=0 

(4) 

with q = 2tt k/ N, k = 0, 1, . . . , N/ 2 and C q being a 
yet arbitrary constant. Thus, a system is said to behave 
diffusively at some wave number q, i.e., on some length 
scale I = 2ir/q, if the corresponding modes F q relax 
exponentially. 

For our purposes, the comparison of the resulting 
quantum dynamics with Eq. (4), it is convenient to 
express the modes F q as expectation values of mode 
operators <!>,, 
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N-l 

F q (t) = Tr{ p(t) $ 9 } , *, = C q J2 cos (M n M , 

(5) 

where the C q are now chosen such that Tr{ } = 1. 

Our approach for the analysis of the F q is based on 
the TCL projection operator technique, see Ref. [8]. 
For the application of this technique we have to define 
a suitable projection superoperator P. Here, we choose 

Pp(t) = T*{p(t)*,}*,=F,(t)*,. (6) 

For initial states p(0) which satisfy P p(0) = p(0) 
[which essentially means that we consider the decay of 
harmonic density waves] the method eventually leads 
to a differential equation of the form 

F q (t) = (/3 2 r a ,,(t) + /3 4 r 4 ,,(t) + . . .) F q (t) . (7) 

V v ' 

r,(t) 

Note that p(0) is not restricted to any energy sub- 
spaces and accordingly corresponds to a state of high 
temperature. Apparently, the dynamics of F q is con- 
trolled by a time-dependent decay rate T q (t). This rate 
is given in terms of a systematic perturbation expan- 
sion in powers of the inter-layer coupling. (Concretely 
calculating the Fi i9 reveals that all odd orders vanish 
for this model.) At first we concentrate on the analy- 
sis of Eq. (7) to lowest (second) order, however, below 
considering the fourth order will account for localiza- 
tion. The TCL formalism yields T2,q(t) = /„ dr f q (r), 
where f q (r) denotes the two-point correlation function 

f q (r) = P 2 Tr{ [ V(t), F q ] [ V(t'), F g ]},T = t-t'. 

(8) 

Here and in the following the time-dependencies of op- 
erators are to be understood w.r.t. to the Dirac pic- 
ture. A rather lengthy but straightforward analysis 
shows that Eq. (8) significantly simplifies under the 
following assumption: The autocorrelation functions 
Tr{ Vp{t) Vft(i') } of the local interactions w M should 
only depend negligibly on the layer number fj, (dur- 
ing some relevant time scale). Simple numerics indicate 
that this assumption is well fulfilled (for the choices of a 
discussed here), once the layer sizes exceed ca. 30 x 30. 
Hence, first investigations may be based on the consid- 
eration of an arbitrarily chosen junction of two layers, 
the interaction in between we label by fj, = 0, i.e., we 
may consider vo in the following. Exploiting this as- 
sumption reduces Eq. (8) to 



/,(t)«-W7(t), (9) 
W = 2 (1 - cos q) 2 , /(t) = l/n Tr{ v (t) v (t') } . 

(Note that the above approximation is exact for iden- 
tical subunits, see Ref. [9].) 

Direct numerical computation shows that /(r) looks 
like a standard correlation function, i.e., it decays com- 
pletely before some time tc- Of primary interest surely 
is T2(t) = / 'dr/(r). Numerics indicate that neither 
to nor 7 = f^t), t > tc [the area under the initial 
peak of /(t)] depend substantially on n (again for n > 
30). Thus, both 7 and tc are essentially functions of a. 
According to all the above findings, an approximative 
evaluation of Eq. (7) to second order reads 

F q (t) » 2 r 2 ,,(t) F q (t) , p 2 r 2 ,,(t) w -WY 2 (t) . 

(10) 

This implies for t > tc 

F q (t)*-W>yF q (t), t r = 1/(Wj). (11) 

The comparison of Eq. (11) with (4) clearly indicates 
diffusive behavior with a diffusion constant k = /3 7. 
Due to the independence of 7 from n, N the pertinent 
diffusion constant for arbitrarily large systems may be 
quantitatively inferred from the diagonalization of a 
finite, e.g., "30 x 30" -layer. 
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Fig. 2. Time evolution of modes F q with q = 2n/N, the 
longest wavelength. Parameters: n = 30, a = 1, N = 10, 
@ = 0.24 (Inset: N = 42, /3 = 1). In both cases the numer- 
ical result (crosses) is an exponential decay which clearly 
indicates diffusive transport behavior and well agrees with 
the TCL2 result (continuous curve), see text for details. 

To check this theory, we exemplarily present some 
results here. For, e.g., a = 1 and n = 30 we numerically 
find tc w 10 and 7 w 2.9. Thus, additionally choosing 
j9 = 0.24 and considering the longest wavelength mode 
in a iV = 10 system (q — 7r/5), we find Wj w 0.064 
[cf. Eq. (11)]. This corresponds to a ratio tr/ tc ~ 1.6, 
that is, tr > tc, which justifies the replacement of 
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Fig. 3. Time evolution of a mode F q with q = 2tt/N, the 
longest wavelength (Inset: q = n, the shortest wavelength). 
Parameters: n = 30, a = 1, AT = 10, /3 = 0.08. The TCL2 
result (continuous curve) fails to describe the dynamics 
correctly for sufficiently long wavelengths. 

Eq. (10) by (11) [see also the discussion of this issue in 
the following paragraph] . And indeed, for the dynam- 
ics of F q (t) we get an excellent agreement of the theo- 
retical prediction based on Eq. (11) with the numerical 
solution of the full time-dependent Schrodinger equa- 
tion (see Fig. 2). Note that this solution is obtained by 
the use of exact diagonalization. Naturally interesting 
is the "isotropic" case of f3 — 1. Keeping a = 1, one 
has to go to the longest wavelength in a N = 42 system 
in order to keep the W from the former example un- 
changed. If our theory applies, the decay curve should 
be the same, which indeed turns out to hold (see inset 
of Fig. 2) . Note that the integration in this case already 
requires approximative numerical integrators like, e.g., 
Suzuki- Trotter decompositions [13]. A numerical inte- 
gration of systems with larger N rapidly becomes un- 
feasible but an analysis based on Eq. (11) may always 
be performed. 

So far, we characterized the dynamics of the diffu- 
sive regime. We turn towards an investigation of its size 
now. Obviously, the replacement of Eq. (10) by (11) is 
only self-consistent for tr > tc, i.e, if the relaxation 
time is larger than the correlation time. This will possi- 
bly break down for some large enough q (small enough 
I), which then indicates the transition to the ballistic 
regime. Since the crossover is expected at tr « tc, 
we may hence estimate the maximum diffusive <j max as 
[cf. Eq. (11)] 

Wmax = VF(g max ) = l/(rc 7) • (12) 

It turns that in the regime of W ^ ITmax a description 
according to Eq. (10) still holds. However, in this case 
T2 (t) is no longer essentially constant but linearly in- 
creasing during the relaxation period. This corresponds 



to a diffusion coefficient k which increases linearly in 
time, which in turn is a strong hint for ballistic trans- 
port (cf. Ref. [9]). Thus, this transition may routinely 
be interpreted as the transition towards ballistic dy- 
namics which is expected on a length scale below some 
mean free path. 

In the following we intend to show that, in the limit 
of long wavelengths, it is the influence of higher order 
terms in Eq. (7) that describes the deviation from dif- 
fusive dynamics. To those ends we consider L, the ratio 
of second order to fourth order terms 



L(t) 



(13) 



/3 2 r 2 , q (t) • 

Whenever L(t) <C 1, the decay is dominated by the 
second order T2, q (t), which implies diffusive dynam- 
ics. It turns unfortunately out that the direct numeri- 
cal evaluation of Y4,, q {t) is rather involved. However, a 
somewhat lengthy calculation based on the techniques 
described in Ref. [14] shows that, for small r 4i9 (£), the 
fourth order term assumes the same scaling in (3, q as 
the second order term and may be approximated as 

/3 4 r4, q (t)~w 2 f 4 (t), (14) 

i 

f 4 (t) = t [l/n 2 ]T(! dT<i|«o(Q«o(OW)-f 2 (t) 2 ] , 
i o 

where \i) are eigenstates of Ho, i.e., F 4 (i) may be eval- 
uated from considering some "representative junction" 
of only two layers, just as done for 1*2 (t). The calcula- 
tion is based on the fact that the interaction features 
Van Hove structure, that is, V 2 essentially is diago- 
nal. We intend to give the details of this calculation 
in a forthcoming publication. Here, we want to con- 
centrate on its results and consequences. [We should 
note that all our data available from exact diagonaliza- 
tion is in accord with a description based on Eqs. (7), 
(10) and especially (14). We should furthermore note 
that r 4 (t), other than 1*2 (i), scales significantly with 
n, which eventually gives rise to the n-dependence in 
Fig. 4.] With Eq. (14) we may rewrite Eq. (13) as 



L(W,t) = W : 



f 4 (t) 
f 2 (t) 



(15) 



This ratio turns out to be a monotonously increasing 
function in t, which is not surprising, since lower order 
terms in Eq. (7) are expected to dominate at shorter 
times. Thus, no visible deviation from the (diffusive) 
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Fig. 4. Numerical results for the measure \ with respect to 
the amount of disorder a and the inverse layer size 1/n 2 . 
An absolute minimum Xmin ~ 0.02 is found at a Ri 0.5 in 
the limit of n — > oo. Note that only 10% of the whole area 
is extrapolated (the area in front of the thick line). 

second order description arises, as long as the decay is 
"over", before L(t) reaches some value on the order of 
a fraction of one. Since the decay time scale is given by 
tr, we are interested in L(W, tr). If L{tr) is on the or- 
der of one, the dynamics of the corresponding density 
wave in the corresponding model must exhibit signifi- 
cant deviations from diffusive, exponential decay. Be- 
cause (apart from the model parameters n, a) tr only 
depends on W [cf. Eq. (11)], we may now, exploiting 
Eq. (15), reformulate L(W,tr) only as a function of 
W and the model parameters n, a but without any ex- 
plicit dependence on f3, q: 

L[W,tr(W)]=R(W) (16) 

We call the above reformulation R(W). It turns out 
that R(W) decreases monotonously with W such that 
the minimum W for which R(W) < 1 holds may be 
found from J?(W m in) = 1. This Wmin corresponds to 
the maximum wavelength beyond which no diffusive 
behavior can be expected. Due to the fact that Ti(t), 
T2(t) and tr(W), tq are numerically accessible, W m in, 
W mK can be computed for a wide range of model pa- 
rameters n, a. In Fig. 4 we display the ratio \ = 
Wmin/Wmax as a function of those model parameters. 
With the approximation W « /3 2 q 2 this ratio allows 
for the following interpretation: 

Wmin Qmin ^min f -\ n\ 

X = W = 12— 

•'max ymax l max 

Hence ^fx (which no longer depends on 0) may be 
viewed as the ratio of the shortest to the longest dif- 
fusive wavelength, the smaller it is, the larger is the 
diffusive regime. 



Obviously, for each layer size n there is some disor- 
der that "optimizes" the diffusive regime (minimizes 
X). But, however, for n = 30 (back of Fig. 4) we find 
y/Xmin ~ 1/3 at this optimum disorder, which indi- 
cates about one diffusive wavelength. Exactly those re- 
spective wavelengths have been chosen for the exam- 
ples in Fig. 2 and the inset in Fig. 3, but not in Fig. 3 
itself. For all a and up to n — 100 (which is about 
the limit for our simple numerics) x clearly appears 
to be of the form x(cr,n) = A(a)/n 2 + B{a). Extrap- 
olating this l/n 2 -behavior yields a suggestion for the 
infinite model n — oo (front of Fig. 4). According to 
this suggestion, we find yxmin ~ 1/7, again at opti- 
mum disorder. This indicates a rather small regime of 
diffusive wavelengths, even for the infinite system. We 
would like to repeat that these findings apply at infi- 
nite temperature, i.e., the above small diffusive regime 
is characterized by the fact that the dynamics within 
it are diffusive at all energies. 
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